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(54) Method and system for generating a diagnostic index 

(57) A method of generating a diagnostic index for 
quantifying the presence or absence of a phenomena 
which is represented by spectral values generated from 
acquired electrical signals of the phenomena, compris- 
es the steps of: 

sorting said spectral values into predetermined bins 
of ranges of spectral values; 
summing all of the spectral values in each bin; 
multiplying the sum of spectral values in each bin 
by a predetermined coefficient to obtain a bin prod- 
uct; 

summing said bin products to obtain a diagnostic 
index which represents a degree of presence or ab- 
sence of said phenomena. 
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Description 

BACKGROUND OF THE INVENTION 

5 [0001] The present invention relates to a method and system for generating a diagnostic index for quantifying the 
presence or absence of a phenomena which is represented by spectral values generated from acquired electrical 
signals of the phenomena. For example, such a method may be used in connection with a computer-based biopotential 
diagnostic system and method for quantitatively determining, in a noninvasive manner, cerebral phenomena that can 
be ascertained by analyzing the properties of cerebral electrical activity. 

10 [0002] Despite a considerable expenditure of time and effort, current approaches to the quantitative, noninvasive 
assessment of cerebral electrical activity, as displayed in an electroencephalographic (EEG) waveform, have not been 
successful in fully extracting all of the information which is present in this complex waveform. A great need remains 
for an accurate, sensitive, reliable, and practical neurological profiling technology. In particular, contemporary intraop- 
erative EEG monitoring techniques have not been widely adopted due to their inherent limitations. Indeed, large num- 

15 bers of medical malpractice suits are believed to be related to post-anesthesia morbidity and mortality, and if such EEG 
monitoring techniques were reliable they certainly would have been adopted. 

[0003] A number of devices known in the prior art are capable of tracking cerebral activity qualitatively. Techniques 
involving the use of the classical conventional analog EEG are restricted to analyses in the time domain, and require 
considerable training for adequate interpretation. Moreover, since the resolution of the human eye at standard EEG 
20 tracing speeds is limited, much of the fine structure of the EEG is invisible. Thus, visual EEG assessment is better 
characterized as an art rather than a science. 

[0004] The use of frequency (power spectrum) analysis of the EEG in the 1 960s introduced the notion of some basic 
processing of the signal prior to visual inspection and led to the application of frequency analysis of the EEG to various 
cerebral monitoring problems. In the past 25 years, over 100 papers have been published in the medical literature 

25 describing applications of power spectral analysis for purposes such as assessing the depth of anesthesia and cerebral 
ischemia under various intraoperative conditions. United States Patent No. 4,557,270 issued to John also describes 
the use of power spectrum analysis to evaluate cerebral perfusion during open heart surgery. Several recent studies, 
however, have shown many deficiencies in the use of power spectral analysis to monitor cerebral perfusion and to 
determine postoperative neurological outcome. In addition, neither power spectrum analysis nor any other monitoring 

30 technique has been shown to be reliable, demonstrated by the fact that the Harvard Medical School Anesthesia Mon- 
itoring Standard does not include any type of intraoperative neurological monitoring, due, in all likelihood, to the com- 
plexity of interpreting raw EEG data and the unreliability of existing automated systems utilizing power spectrum or 
time-domain analytic techniques. 

[0005] The discharge of thousands of bioelectricaily active cells in the brain, organized in larger, interacting neural 

35 centers contributes to the formation of an electrical signal with a wide frequency spectrum that is rich in harmonics and 
extremely complex dynamics. Embedded in that signal is information regarding frequency content, nonlinearities, and 
phase relationships arising from the complex neuronal firing patterns that take place. Such firing patterns change 
constantly making the stastisticai properties of the EEG signal highly nonstationary, Because of the complexity of the 
EEG signal, conventional time and frequency modes of analysis have not been able to fully profile its behavior. This 

40 may be one of the reasons for the limited success of such approaches. 

[0006] In the Fourier transform of the second order autocorrelation function (the power spectrum), processes are 
represented as a linear summation of statistically-uncorrelated sine-shaped wave components. Contemporary ap- 
proaches to monitoring the EEG by means of the power spectrum have thus suppressed information regarding non- 
linearities and inter-frequency phase relationships and are of limited utility in representing the EEG's dynamic structure. 

45 [0007] Because the EEG is highly dynamic and nonlinear, the phase relationships within the EEG are the elements 
most likely to carry diagnostic information regarding cerebral function. The Fourier transform of the third order auto- 
correlation function, or autobispectrum, is an analytic process that quantifies deviation from normality, quadratic non- 
linearities and inter-frequency phase relationships within a signal. The Fourier transform of the third order cross-cor- 
relation function, or cross bispectrum, is an analytic process that provides similar information for two signals. We can 

50 generalize these techniques by defining the Fourier transform of the nth-order auto/cross correlation function, or the 
n-1 order auto/cross spectrum, as an analytic process that contains information regarding deviation from normality, as 
well as n-1 order nonlinearities and inter-frequency phase relationships in a signal. Auto/cross spectra beyond the 
bispectrum will be referred to as higher-order spectra. 

[0008] Autobispectrum analysis techniques have been applied to the EEG signal to demonstrate the basic bispectral 
55 properties of the conventional EEG. Such studies have also been conducted to search for differences between the 
waking and sleeping states. Autobispectrum analysts and power spectrum analysis have also been used in an attempt 
to show that the EEGs of monozygotic twins are similar in structure. United States Patents No. 4,907, 597 and 5,01 0,891 
issued to Chamoun describe the use of auto/cross bispectrum analysis of the EEG to evaluate cerebral phenomena 
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such as quantifying depth and adequacy of anesthesia, pain responses induced by surgical stress, cerebral ischemia, 
consciousness, degrees of intoxication, ongoing cognitive processes and inter-hemispheric dynamic phase relations. 
[0009] To date, no one has used auto higher-order spectrum or cross higher-order spectrum analysis for neurological 
diagnoses or monitoring of the cerebral phenomena described above. 

5 [001 0] A common problem in analyzing the data generated by any of the spectral techniques discussed above is the 
fact that the EEG's frequency distribution may dramatically change under relatively stable physiological conditions. 
Such changes will lead to changes in the power spectrum, bispectrum, and higher order spectra at the corresponding 
frequencies. For example, when hypnotic anesthetic agents are administered in low to medium concentrations, there 
is a substantial increase in the EEG activity in the 12-18 Hz frequency band. High doses of the same agents will lead 

10 to a sudden reduction in activity in the 1 2-18 Hz band and increase in activity in the 0.5-3.5 Hz band, followed by burst 
suppression at extremely high concentrations. A frequency-based analysis that uses the 12-18 Hz frequency band to 
track the patient's anesthetic depth during the administration of a hypnotic agent will provide a misleading assessment 
of the patient's depth when the shift in activity from high to low frequency occurs. Such transitions are even more 
complicated when a mixture of anesthetic agents is used. 

w [0011] It is an object of the present invention to provide an effective method and system of generating a diagnostic 
index from spectral values generated from acquired electrical signals. 

[0012] According to the invention this is achieved by a method with the features of claim 1 and by a system with the 
features of claim 2. 

[001 3] These and other objects and features of the present invention are more fully explained by the following detailed 
20 description and figures. 

BRIEF DESCRIPTION OF THE FIGURES 

[0014] 

25 

Fig. 1 is a schematic view of a system for detecting cerebral phenomena in a non-invasive manner; 

Fig. 2 is a schematic view of a 19 channel EEG data acquisition and analysis system utilized in the system of Fig. 1 ; 

30 Fig. 3 is a schematic view of the microcomputer used to display the EEG power spectrum and bispectrum higher- 

order spectrum in the system of Fig. 1 ; 

Fig. 4 is a schematic view of the processing operations performed by the system of Fig. 1 ; 

35 Fig. 5 is a flow chart of the operations of the monitor module shown in Fig. 4; 

Figs. 6(a) - 6(c) are views of sample display representations of diagnostic index generated by the system of Fig. 1 ; 

Fig. 7 is a flow chart of the operations of the acquisition and EEG raw data management module of the system 
40 shown in Fig. 4; 

Fig. 8 is a flow chart of the frequency-domain-based method for producing autobispectrum, cross bispectrum, auto 
power spectrum, or cross power spectrum used by the system of Fig. 1 ; 

45 Fig. 9 is a flow chart of the parametric based method for producing autobispectrum, cross bispectrum, auto power 

spectrum, or cross power spectrum in the system of Fig. 1 ; 

Fig. 10(a) is a graph showing a bispectral density array generated by the system of Fig. 1; 
so Fig. 10(b) is a graph showing a biphase array generated by the system of Fig. 1; 

Fig. 10(c) is a graph showing a bicoherence array generated by the system of Fig. 1; 

Fig. 10(d) is a graph showing an array of square root of real triple product generated by the system of Fig. 1; 

55 

Fig. 11 is a flow chart of the operations of the diagnostic index derivation module shown in Fig. 4; 
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DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS 

[001 5] Referring to Fig. 1 , the apparatus for noninvasive^ detecting cerebral phenomena includes a 1 9 channel EEG 
data acquisition and analysis system 12 connected to a microcomputer 18. 
s [001 6] The EEG leads are connected to a patient's head 1 4 by a set of surface electrodes 1 3. The international 1 0/20 
electrode system is preferred. The EEG signals are detected by the electrodes and transmitted over a patient cable 
16 to the EEG data acquisition and analysis system 12. 

[0017] The data acquisition and analysis system 1 2 filters and amplifies the EEG waveforms. Commonly used digital 
signal processing (DSP) techniques are applied to digitize, to low-pass filter (50 Hz), and to decimate the signals. 

10 Power spectral, bispectral, and higher-order spectral processing can then be performed. 

[0018] The system 12 generates all power spectrum, bispectrum, and higher-order spectrum arrays. These arrays 
are then used in conjunction with clinically predetermined coefficient arrays to produce diagnostic indices. These indices 
are sent to the host computer 18 and are displayed on the graphics display 20. Printed output of the diagnostic index 
is also available on the hard copy output device 22 which is connected to the microcomputer 18. The operator interacts 

is with the acquisition and analysis components of the system by means of a user input device 24 with feedback on the 
graphics display 20. 

[0019] The 19 channel data acquisition and analysis system 12 is shown in greater detail in Fig. 2. The EEG surface 
potential, detected by surface electrodes 13 mounted on the patient's head 14, passes through an electros urgery 
protection circuit 30, a defibrillator protection circuit 32, and an amplifier/filter circuit 36 before being passed on to the 

20 multi-channel analog to digital converter 38. 

[0020] The electrosurgery protection circuit 30 includes a radio frequency (rf) filter, which limits the rf current through 
the patient leads 16 to less than 100 microamperes and thus protects the patient 15 from rf burns and protects the 
amplifiers 36 from damage resulting from exceeding the absolute maximum input voltage specified by the manufacturer. 
This circuit can be an LC circuit consisting of a generic inductor connected in series to a generic capacitor which is 

25 then connected to ground. 

[0021] The defibrillator protection circuit 32 limits the voltage to the amplifiers 36 to a safe level when a defibrillator 
is applied to the patient 15 and discharged. This circuit can consist of a generic resistor connected, in series with the 
signal path, to a neon light bulb or other surge suppression device which is then connected to ground. 
[0022] The amplifier/filter circuitry 36 is controlled by the microcomputer 18 for gain and filtering levels which may 

30 be adjusted by the operator. Preferred gain and filtering settings are discussed below. This circuit section consists of 
three stages. The first is a pre-amplifier stage that can be assembled using a wide variety of high-impedance pre- 
amplifiers such as those sold by National Semiconductor, Sunnyvale CA. The second is a stage composed of pro- 
grammable filters which will allow an adjustable band pass cutoff to be selected anywhere in the range of 0.1 Hz to 4 
KHz. The filters can be designed using components from Frequency Devices, Haverhill MA. The third stage is composed 

35 of programmable amplifiers which can be assembled from operational amplifiers used in conjunction with a multiplying 
digital to analog (D/A) converter. Both components can be obtained from National Semiconductor. The multiplying D/ 
A is used to set the gain to the appropriate levels requested by the microcomputer 1 8. 

[0023] The high impedance pre-amplifier of each channel will saturate to either the positive or negative supply voltage 
if the input of the pre-amplifier is not terminated. This will lead to large positive value or a large negative value at the 
40 output of amplifier/filter section 36. Such values will be used to identify lead failure. 

[0024] The output of all 1 9 channels of the amplifier/ filter 36 is fed to the multi-channel A/D converter 38 which is 
controlled by an input processor 44 for sampling rate settings. The analog signals are converted to digital data format 
suitable for input to the input processor 44. A/D converters sold by Analog Devices, Norwood MA can be used for this 
purpose. 

45 [0025] The multi-channel A/D converter 38 is optically coupled to the input processor 44 by optical isolator 40. All 
control lines to the A/D converter 38 are also optically isolated by optical isolator 42. Any optical isolator can be used 
for this purpose. 

[0026] All DC power lines connected to the amplifiers 36 and A/D converter 38 are also isolated from the AC power 
line with a DC/DC convenor 43 in order to provide complete patient isolation from ground. DC/DC converters available 

so from Burr Brown can be used for this purpose. 

[0027] The basic instructions for controlling operation of the input processor 44 are stored in a read only memory 
(ROM) 46. The random access memory (RAM) 48 is used as a buffer memory for data and a portion of the RAM 48 
can also be used as program memory when a control program is being downloaded from the microcomputer 1 8. The 
input processor 44 has a bus 50 to communicate with its RAM 48 and ROM 46 and a separate bus 55 for communicating 

55 with the microcomputer 1 8. 

[0028] The memory architecture of the calculation processor is similar to that of the input processor. The basic in- 
structions for controlling operation of the calculation processor 52 are stored in a read only memory (ROM) 54. The 
random access memory (RAM) 56 is used as a buffer memory for data and a portion of the RAM 56 can also be used 
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as program memory when a control program is being downloaded from the microcomputer 18. The calculation proc- 
essor 52 has a bus 58 to communicate with its RAM 56 and ROM 54 and uses the bus 55 for communicating with the 
microcomputer 18. 

[0029] The A/D converter 38 acquires the data at high speed and filtering is done by the input processor 44 to exclude 
5 frequencies outside the region of interest. The input processor simultaneously decimates the sampling rate of the input 
data to a lower sampling rate. The input processor 44 transfers the filtered and decimated data stream to the micro- 
computer 18 for display of the raw input signals via the data bus 55 and buffers 60 to the microcomputer data bus 40. 
The input processor 44 also transfers the data to the calculation processor 52 for calculation of power spectrum and 
higher-order spectrum characteristics of the input signals via a serial communication interface 51. The calculation 
w processor 52 calculates power spectrum and higher-order spectrum characteristics of the input data and produces 
diagnostic indices from the calculated power spectrum and higher-order spectrum data. The input processor can be 
any general purpose DSP processor such as the ADSP-2101 sold by Analog Devices, Norwood MA. The calculation 
processor is a floating-point DSP processor in the preferred embodiment such as the TMS320C30 sold by Texas 
Instruments, Dallas, TX. 

15 [0030] The host or microcomputer 18 of Fig. 1 is shown in greater detail in Fig. 3. The entire microcomputer system 
runs under control of a microprocessor 62 with the program memory stored in ROM 64. The RAM 66 is used for storage 
of intermediate data. The storage device 84 can be a Winchester disk or a large block of RAM or any other storage 
medium. It is used for storage of clinical information and can be used for archiving patient data. 
[0031 ] The microcomputer 1 8 may contain a math coprocessor 70 which is connected directly to microprocessor 62. 
20 The math coprocessor 70 is used for scalar and graphic calculations. A graphics controller 72 operating under program 
control of the microprocessor 62 drives a graphics display 20. An interface port 74 provides the connection from the 
microcomputer bus 40 to the user interface device 24. The user interface device 24 may be a keyboard, a pointing 
device or a keypad or any combination of these or similar devices. The interface port 74 can also provide a connection 
between the microcomputer and an external evoked potential stimulating device. This connection will allow the micro- 
ns computer to trigger a stimulus or easily identify the onset of an independently triggered stimulus. 

[0032] Operator control of the entire acquisition, analysis and display procedure is controlled by the user interface 
device 24 with feedback on the graphics display 20. The data bus 40 can be used to send control data to the 1 9 channel 
data acquisition system 12 (e.g., filtering, gain, sampling rate, start/stop acquisition, perform self diagnostics) and to 
receive EEG data from the system, as well as to download program data to the system. A serial or parallel port 78 is 
30 provided to drive a hard copy output device 22 for printing desired diagnostic indices. 

[0033] Referring now to Fig. 4, a block diagram of the system operations is described. As mentioned above, the 
system for noninvasiveiy detecting cerebral phenomena computes dynamic phase and density relations of EEG signals 
from a preselected number of leads. Single-valued diagnostic indices are then generated from the data arrays by using 
clinically predetermined coefficient arrays. The results are quantitative indices useful for analyzing cerebral electrical 
35 activity as it relates to, for example, the assessment of depth and adequacy of anesthesia, pain responses during 
surgical stress, cerebral ischemia, cerebral hypoxia, level of consciousness, degree of cerebral intoxication, altered 
evoked potential responses, and normal or abnormal cognitive processes that include but are not limited to Alzheimer's 
disease and HIV-related dementias. 

[0034] The monitor module 402 handles the overall operations of the system via integration of data and process 
40 information from the user interface module 404, acquisition and raw EEG data management module 406, power spec- 
tral, bispectral and higher-order spectral processing module 408, and the diagnostic index derivation module 410. A 
detailed illustration of module 402 can be found in Fig. 5. 

[0035] The operator controls and interacts with the system during the course of a procedure through the user interface 
and display management module 404. This interaction includes, but is not limited to, entry of information regarding the 
45 patient and type of diagnostic procedure underway; lead and acquisition settings; continuous display of acquisition 
status, lead integrity, and diagnostic indices corresponding to regions probed by each electrode; and requests for 
printing and archiving results to the storage device. Module 404 directly interacts with the monitor module 402. The 
operations handled by module 404 can be achieved under a commercially available environment such as Microsoft 
Windows. 

50 [0036] The acquisition and raw EEG data management module 406, handles all of the raw EEG data checking and 
processing prior to power spectrum, bispectrum, and higher-order spectrum analysis. This includes, but is not limited 
to, continuous acquisition of EEG data and the verification of its integrity; preparation of all unipolar EEG data for auto/ 
cross power spectral, bispectral, and higher-order spectral processing. Module 406 directly interacts with the monitor 
module 402. A more detailed description of module 406 is provided below in connection with Fig. 7. 

55 [0037] The power spectral, bispectral, and higher-order spectral processing module 408 controls the generation of 
all data arrays for power distribution, dynamic phase relations, and power coupling within the EEG. This information 
can be obtained by computing the auto/cross power spectrum, bispectrum, and higher-order spectra using either an 
FFT-based or parametric-based approach. The tasks performed by this module include, but are not limited to: Fourier 
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transformation and the generation of power spectra; auto/cross bispectral density and higher order density generation; 
auto/cross bicoherence and higher order coherence generation; auto/cross bispectral real product and higher-order 
real product generation; and auto/cross biphase and higher-order phase generation. Module 408 directly interacts with 
the monitor module 402. A more detailed description of module 408 is provided below in connection with Figs. 8 and 9. 

5 [0038] The diagnostic index derivation module 41 0 generates the data values used in the diagnostic process. The 
task includes, but is not limited to, sorting the values in the frequency band of interest for each of the required power 
spectrum, bispectrum, or higher-order spectrum arrays; dividing each of the sorted arrays into bins representing por- 
tions of the distribution histogram of the sorted data (i.e. top 0-5%, top 5-10% as well as bottom 5%, etc.); summing 
the values in each bin to create a single number variable; creating a diagnostic index by multiplying the resultant sorted 

10 values from auto/cross power spectrum, bispectrum, and higher-order spectrum arrays by clinically predetermined 
coefficients; and summing all variables that have been multiplied by a coefficient to create a final univariate diagnostic 
• index. Module 41 0 directly interacts with the monitor module 402. A more detailed description of module 41 0 is provided 
below in connection with Fig. 11 . 

[0039] A schematic of the operation of the monitor module 402 is shown in Fig. 5. In initializing step 502, the data 
15 arrays are filled with the most recent 63 seconds of raw digitized EEG signals, and the power spectrum, bispectrum, 
and higher-order spectrum data for each lead are initialized to zero. The data files required for storage and files con- 
taining data bases required for the computation of diagnostic indices, are also opened in the initializing step 502. 
[0040] In step 504 the system requests the information required to start the acquisition and diagnostic process from 
the user via the user interface module 404. This requested information includes patient descriptive statistics (sex, age, 
20 clinical symptoms, etc.), type of diagnostic procedure to be conducted, the leads used for auto power spectrum, bis- 
pectrum, and higher-order spectrum analysis as well as the leads to be used for cross power spectrum, bispectrum, 
and higher-order spectrum analysis. 

[0041] In its default mode of operation the system continuously monitors the depth and adequacy of anesthesia and 
any pain responses during surgical stress using a default autobispectrum database. Default band-pass filtering is 
25 performed, passing the range 0.5 to 50 Hz; the default sampling rate is set at 1 28 samples per second; and the default 
gain is set at 5000 for each lead. The following discussion and description of the preferred embodiments will emphasize 
autobispectral processing performed on EEGs from specific electrode sites that best provide depth of anesthesia in- 
formation. Other modes of operation will be described more generally. 

[0042] According to the international 10/20 electrode system, the 19 EEG signals that can be acquired using the 
30 system are: Fp1 , Fp2, F7, F3, Fz, F4, F8, T3, C3, Cz, C4. T4, T5, P3, Pz. P4, T6, 01 , and 02 (A1 or A2 for reference). 
[0043] In order to perform auto power spectrum, bispectrum, and higher-order spectral analysis, one signal is re- 
quired. This signal can be measured directly from any of the above electrodes or it can be synthesized by linearly 
combining signals from two or more EEG leads. For example, two analog signals can be subtracted from each other 
using a differential amplifier to yield a third signal. The same operation can be performed on the two digitized signals 
35 using numerical subtraction. The auto power spectrum data will provide information regarding the power distribution 
within the signal; the autobispectrum data will provide information regarding deviation from normality, quadratic non- 
linearities and inter-frequency phase relationships within the signal; finally, auto higher-order spectrum data will provide 
information regarding deviation from normality, higher-order nonlinearities, and inter-frequency phase relationships 
within the signal. Such processing will determine if the signal is made up of independent wave components or whether 
40 certain frequencies are simply harmonics of nonlinearty interacting fundamentals. Cerebral phenomena that alter the 
nonlinear frequency structure of the signal at the location probed by the electrode are best quantified by autobispectrum 
and higher-order spectrum type approaches. 

[0044] In order to perform cross power spectrum, bispectrum, and higher-order spectrum analysis, two signals are 
required. The two signals can be measured directly from any of the above electrodes or either of the two signals can 

45 be synthesized by linearly combining two or more of the EEG leads as described earlier. The cross power spectrum 
data will provide information regarding the power correlation between the two signals. The cross bispectrum data will 
provide information regarding deviation from normality, quadratic nonlinearities, and inter-frequency phase relation- 
ships between the two signals. Finally, cross higher-order spectrum data will provide information regarding deviation 
from normality, higher-order nonlinearities, and inter-frequency phase relationships between the two signals. Such 

50 processing will determine if the frequencies in signal "X" are independent or whether they are harmonics of fundamen- 
tals present in signal °Y°. This provides a better characterization of the relationship between two signals originating 
from separate regions of the cortex. Cerebral phenomena that alter nonlinear frequency relations between the various 
regions of the cortex are best quantified by cross bispectrum and cross higher-order spectrum approaches. 
[0045] Since the effects of anesthesia are reflected by more global changes in the EEG, the preferred embodiment 

55 will use six signals to illustrate the operation of the system using autobispectrum analysis for the monitoring of the 
depth of anesthesia. The six signals are derived from the following electrode placements: left and right frontal (FL7FR) 
signals are derived from (Fp1-Cz) and (Fp2-Cz) respectively; left and right parietal (PL/PR) signals are derived from 
(P3-Cz) and (P4-Cz) respectively; left and right fronto-parietal (FPIVFPR) signals are derived from (Fp1-P3) and 
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(Fp2-P4) respectively. 

[0046] In step 506, a new one-second buffer of raw EEG data is acquired. The system performs artifact detection 
on the new one-second buffer and properly updates all data arrays. Any transmission of artifactual data is displayed 
to the operator in order to invoke the operator into correcting the problem. 
5 [0047] The system, in step 508, computes auto power spectrum and autobispectrum arrays for the signals FL, FR, 
PL, PR, FPL, FPR. Other signals may, of course, be used for auto/cross power spectral, bispectral, and higher-order 
spectral processing. Two different approaches for power spectrum, bispectrum, and higher-order spectrum computation 
will be discussed below with reference to Figs. 8 and 9. 

[0048] In step 510, the single-valued diagnostic indices from all generated auto/cross power spectrum, bispectrum, 
10 and higher-order spectrum arrays are computed. The clinically predetermined coefficient arrays for the auto/cross 
power spectrum, bispectrum, and higher-order spectrum arrays are used for the diagnostic index computations. The 
generation of the coefficient arrays is discussed later. The system instantaneously displays, in step 512, all computed 
diagnostic indices for all signals being analyzed. In step 514, the system checks for an exit request, and if such a 
request has not been made, the system repeats steps 508 through 514. In step 516, requested printouts are produced, 
15 results are stored to a storage device for archival purposes and all files are closed. In step 51 8, the process is terminated. 
[0049] A sample condensed display representation generated by the system is shown in Figs. 6(a) - 6(c). Represen- 
tations of the patient's head are shown on the graphics display in Fig 6(a) and Fig. 6(b). The first illustration Fig. 6(a) 
is divided into nineteen sections each representing the region probed by an electrode. The second illustration Fig. 6(b) 
is divided into three horizontal sections representing combined left and right hemisphere activity probed by a group of 
20 electrodes in that region. The virtual head displayed on the screen may be partitioned as required for a particular 
diagnostic or monitoring application. For example, if a global effect like depth of anesthesia is being tracked, then one 
unified index along with its trend may occupy the whole display area. 

[0050] For head representation Fig. 6(a), each section contains the instantaneous value of the index 602 using EEG 
data acquired from the electrode in that regionfFor head representation Fig. 6(b), each section contains the instanta- 
25 neous value of the computed index 604 using EEG data acquired from several electrodes in that region. Next to each 
index value, a color-coded arrow is used to show the instantaneous change in the direction of the index. The arrow 
will be green if the index is within acceptable limits set as by the operator. The arrow will change to yellow if the index 
moves into a warning zone. A flashing red bar will replace the arrow if the index has a value that is outside the acceptable 
limits set for the patient. 

30 [0051] At the request of the operator, the instantaneous value of the index and its trend for any section can be 
displayed as an enlarged view 606 for closer examination as shown in Fig. 6(c). This will facilitate the examination of 
the patient's status at a distance. Each section will be covered by a large "X" 608 if a lead fails or artifact was detected, 
for any of the leads contributing to the data required to generate the diagnostic index for that region. 
[0052] Referring to Fig. 7, the acquisition and raw EEG data management module 406 will now be described in 
35 greater detail. In step 702, the system checks whether new data is being acquired for the first time. If it is, the acquisition 
system 12 in step 704 is supplied with requested filtering, gain, sampling rate, and lead selection information. The 
default settings are band pass 0.5 - 50 Hz for filtering, 5000 for gain, 128 samples/sec for sampling rate and signals 
from the lead combinations FL, FR, PL, PR, FPL and FPR are acquired. The above settings are quite different when 
the system is analyzing evoked EEG responses rather than continuous EEG signals. Common gain and filter settings 
40 to acquire signals for the various EEG evoked potentials are described below. 

? [0053] EEG evoked potentials are a means by which the sensory areas of the brain and of the central nervous system 
/ may be assayed by detecting responses in the EEG to sensory stimuli. There are three common methods: Pattern- 
/ shift visual evoked potentials (PSVEP) involve a visual pattern that is shown to the patient and changed. For example, 
a strobe light may be flashed or a black and white checkerboard may be reversed (black for white and vice versa). 
45 Brainstem auditory evoked potentials (BAEP) uses a controlled auditory stimulus such as a click produced by a signal 
A . generator. Finally, somatosensory evoked potentials (SEP) employs either physiologic (touch or muscle stretch) or 

electrical stimuli. In all evoked potential methods, electrodes are placed near the appropriate centers of the brain (i.e. 
over the visual cortexes in the case of visual evoked potentials) and EEGs are recorded for a certain period of time 
beginning with the administration of the stimulus. The stimulus is repeated many times and the resulting recordings 
so are averaged (traditionally, in the time domain) so as to eliminate all parts of the EEG signal except that due to the 
v stimulus. In the present invention, a series of power spectrum, bispectrum, or higher-order spectrum arrays, as pro- 
duced from the EEG of the evoked responses, is averaged. 
[0054] For each evoked potential method, different filter and gain settings are used. For example, a range of common 
gain settings for pattern-shift visual evoked potentials is 20,000 to 100,000. A range of common filter settings for 
55 PSVEPs is 1 to 3 Hz for the low end of the band pass and 100 Hz to 300 Hz for the high end. The methods and use 
of evoked potentials are described more fully in Evoked Potentials In Clinical Medicine, by Chiappa 1983. 
[0055] In step 706, the acquisition system 12 acquires one second's worth of new data for all requested leads. Or 
the signal from one complete evoked potential response is acquired if the system is analyzing evoked potentials. The 
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system detects lead failures during the acquisition cycle in step 708 by checking for very large positive or negative 
values. Also in step 708, publicly available algorithms are used to check for artifacts in each lead. In step 710, leads 
that have failed and those producing artifactual data are marked for the monitor module 402. 
[0056] In step 712, the most recent 4-second record for each of the signals is assigned to Xj(t), where Xj(t) is the 

5 individual time series records provided for auto power spectral, autobispectral, and auto higher-order spectral process- 
ing (herein, the time series X s (t) (for all t, for one specific i) is referred to as a record). In situations where cross power 
spectral, bispectral, and higher-order spectral processing is required, the most recent 4-second record from the second 
signal is assigned to Y ( (t). In the preferred embodiment, Yj(t) is set to equal Xj(t) in all cases, since only auto power 
spectrum, auto bispectrum, and auto higher-order spectrum computations are to be performed. The index i denotes 

10 the record number from 1 to 60. If evoked potentials are being analyzed, the most recent complete evoked potential 
response from each signai is assigned to the appropriate X,(t) and Y,(t) as described above. Using evoked potential 
responses as individual records will allow us to average a large number of them in the power spectrum, bispectrum, 
and higher-order spectrum domains. 

[0057] In step 714, a circular buffer mechanism is used for storing the raw EEG for each lead, as well as the auto/cross 
15 power spectrum, bispectrum, and higher-order spectrum arrays for the sixty most recent 4-second Xj(t) and Yj(t) records 
for each lead. The buffer is updated by storing the most recently acquired and processed data in the location of the 
oldest data. Operation of the system returns to the monitor module 402 in step 716. 

[0058] Referring now to Fig. 8, the frequency-domain-based procedures for producing the auto power spectrum, 
autobispectrum, cross power spectrum, or the cross bispectrum will now be discussed. 
20 [0059] In step 802, the system checks whether the computation to be performed requires one signal or two signals. 
Typically, one time series is required to perform autospectrum analysis and two time series are required to perform 
cross spectrum analysis. 

[0060] In step 804, the system sets time records in the following manner in order to proceed with auto power spectral 
or autobispectral computation of the unipolar lead. As these computations require only one signal, the second set of 
25 records (Yj(t)) is set to equal the first set (Xj(t)). As a consequence, the corresponding Fourier transforms of Xj(t) and 
Yj(t), respectively Xj(f) and Yj(f), are also equal: 

Xfi) = Yft)~> XjtfHYjtf) 

30 

where i denotes the record number which, in this embodiment, ranges from 1 to 60. 

[0061] In step 806, time records are set for cross power spectral and cross bispectral analysis using two separate 
time series signals. As a consequence, the corresponding Fourier transforms are not equal: 

35 XjW^YjW-^X^^Yjtf) 

where X s (t) and Yj(t) represent individually derived time series records from two separate regions probed by two or 
more electrodes. 

40 [0062] The fast Fourier transform (FFT) Xj(f) and Y,(f) of each of the 60 selected records for that signai, is computed 
using a standard IEEE library routine (or any other publicly available routine) in step 808. If requested, the series of 
transformed records, Xj(f) and Yj(f), may be each normalized by dividing the value at each frequency by the constants 
C xi and C yi , respectively. These constants are derived separately for each record and each series (either X or Y). The 
constant could be the total power, the largest peak in the spectrum of interest, or some other derivative of X|(f), Xj(t), 

45 Y ( (f), and Yj(t). 

[0063] In step 81 0, the system checks whether the computation to be performed is a power spectrum or bispectrum 
computation. 

[0064] The system computes the auto/cross power spectral density values (PD(f)) in step 81 2 by using the following 
equations where PC(f) is the average complex product for a signai or signal pair: 

50 

M 

PC(fJ = ~ X Xi(f) * ^ (f) 



PD(f) = IPC(f)l 
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where /.(f) is the complex conjugate of Yj(f) (0 < f < N/2) and M is the number of records (60 in the preferred embod- 
iment). The system then returns the requested auto/cross power spectral density array to monitor module 402. 
[0065] If the system is performing a bispectral computation in step 81 4, the system checks whether the computation 
to be performed is an autobispectrum or cross bispectrum computation. 

[0066] Autobispectrum analysis is a special case of crossbispectrum analysis and therefore different rules of sym- 
metry apply. In step 816, the system uses the following equations to determine what ranges of fy and f 2 to use during 
autobispectral computation: 

< 4 S N 

where N equals the number of samples per record (512 = 4 sees per record * 128 samples per sec in a preferred 
embodiment), and 

0 < f 2 < f A 

where and f 2 (also referred to as F-, and F 2 or Frequency 1 and Frequency 2) denote the frequency pairs over which 
bispectrum computation will be carried out. 

[0067] In step 81 8, the following equations are used to determine the range of f, and f 2 for cross bispectrum analysis: 

0 < f, + f 2 Z 5 

f 2 <f, 

where all variables represent the same values as they do for autobispectral analysis, except that for crossbispectral 
analysis X s (f) and Yj(f) represent the Fourier transform of the individually derived time series records from two separate 
regions. 

[0068] In Step 820, the power spectra P xj (f) and P^f) of each of the 60 selected records for that signal are computed 
by squaring the magnitudes of each element of the Fourier transform Xj(f) and Y^f) respectively. 
[0069] The system computes the average complex triple product in step 822 by using the following equations where 
bc^.fg) is an individual complex triple product from one 4-second record and BCtf^) is the average complex triple 
product for all 60 records: 

bc i (f 1 .f 2 ) = X | (f 1 )*X j (f 2 )*v;(f 1 +f 2 ) 
where / / (f 1 +f 2 ) is the complex conjugate of Y i (f 1 +f 2 ), and 

M 

BC(fi,f2) =7, ^bcufi.fz) 
i = 1 

where M is the number of records (60 in the preferred embodiment) 

[0070] The average real triple product is computed in step 824 by using the following equations where br^, f 2 ) is 
an individual real triple product from one 4-second record and BR(f 1t f 2 ) is the average real triple product for all 60 
records: 
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br i (f 1 ,f 2 ) = P xi (f l )*P xi (f 2 )*P yi (f 1 +f 2 ) 



M 

BR(fi / f 2 ) - ^ 5>ri(f:,f 2 ) 

i = 1 

where M is the number of records (60 in the preferred embodiment) 

[0071] In step 826, the array of auto/cross bispectral density values (BD(f 1f f 2 )) are computed using the following 
equation: 

BD(1 v i 2 ) = \BC(1 v 1 2 )\ 

[0072] In step 828, the array of the square roots of the average real triple products (SBR(f 1f f 2 )) are computed using 
the following equation: 

SBR(f 1 ,f 2 ) = [BR(f 1 ,f 2 )] 1/2 

[0073] In step 830, the system computes the array of auto/cross biphase values (cp(f 1 ,f 2 )) using the following equation: 

<p(f n ,f 2 ) = tan' 1 [lm(BC(f 1 ,f 2 ))/Re(BC(f 1 ,f 2 ))] 
0 < <p < 2k (radians) 

[0074] In step 832, the system computes the array of auto/cross bicoherence values (R(f 1 ,f 2 )) using the following 
equation: 

R(f 1 ,f 2 ) = BD(f 1 ,f 2 )/SBR(f 1 ,f 2 ) 



0<R<1 

[0075] In step 834, the system returns the requested auto/cross power spectral density array or auto/cross bispectral 
density, squared-rooted average real triple product, bicoherence, biphase arrays to the monitor module 402. 
[0076] The above frequency-domain-based equations used to compute the auto/cross bispectrum arrays can be 
generalized to compute auto/cross higher-order spectral arrays. This will allow the computation of the trispectrum, 
quadspectrum, etc. Assuming that the arrays for a Kthorder spectrum are to be computed the following equations can 
be used: 

[0077] The average complex Kth order product: 

M 

KC(f:,f 2 fK-l) = ^ E X i (fl) * Xi(f2) * 

i = 1 

* Xi(fK-i) * 3Ti{fl+f2+...+fK-l) 

where M is the number of records (60 in the preferred embodiment) 
[0078] The average real Kth-order product: 
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M 

KR(f:,f 2 f K -l) =5 X p xi(fl) * Pxi(f2) * ... 

* PxivfK-l) * Pyi(fl+f2+. . -+fK-lJ 
w [0079] The auto/cross Klh-order spectral density: 

KD(f v f 2 f K -i) = | KC ( f v f 2-- f K-i)| 

15 [0080] The auto/cross Kth-order coherence: 

RP,, f 2> -..,f K . 1 ) = KD (f v f 2 ff^VIKRtf,, f 2 f^)] 1 * 

0<R<1 

[0081] The auto/cross Kth-order phase: 

25 <p(f 1 ,f 2 , f K .i) ■ tan" 1 [lm(KC(f 1 ,f 2 f K-1 )) 

/Re (KC(f v f 2 f K-1 »] 

0 < <p < 2n (radians) 

[0082] Figure 9 illustrates a parametric-based method for producing the auto power spectrum, autobispectrum, cross 
power spectrum, or cross bispectrum. In steps 902, 904, and 906 the system sets the time series records in the same 
35 manner as described above in steps 802, 804, and 806 respectively. The auto/cross power spectra of Xj(t) and Yj(t) 
are estimated in steps 908, 910, and 912. This estimation method includes two major stages, the autoregressive (AR) 
model order selection and the auto/cross power spectrum computation for Xj(t) and Yj(t). In step 908, the system 
computes two sequences of autocorrelations, {R2 X (m)} and {R 2Y (m)} using the following equation. 

40 

H N-l/nl 
i=l C=0 

45 

z = X or Y, and m = 0, 1 L 



where M is the number of records of each signal (60 in the described embodiment), and N is the number of samples 
so per record (512 in the described embodiment), and L is greater than the largest possible AR filter order (50 in the 
described embodiment). 

[0083] The Final Prediction Errors, FPE x (m) and FPE Y (m) are calculated for all orders, m = 1 , 2, .... L, by performing 
a Levinson recursion function on each autocorrelation sequence in step 910 in order to find the order of the AR filter. 
The order of the AR filters can be determined by finding the location of the minimum of Final Prediction Errors: FPE x (m) 
55 and FPEy(m) respectively, i.e., 

FPE^Q X ) = min [FPE^m)} 
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and 



FPEy{Q Y ) = min {FPE^m)} 

5 

where Q x and Q y are the locations of the minimum values for FPE x (m) and FPE y (m) (respectively) and, consequently, 
the orders of the AR filters of the power spectra Xj(t) and Yj(t) (respectively). 

[0084] Once the orders of the AR filters for auto power spectra are known, the autocorrelation sequences, (R^m)} 
and {R 2 Y( m )}» are entered into a Levinson recursion with order Qx and Qy, respectively, instead of L. The coefficients, 
10 fox. i=0> 1 . — .Qx) and { c iY- i = 0.1 — .Qy}* obtained from the recursion are the coefficients of the AR filters for auto 
power spectra of Xj(t) and Yj(t) respectively. Then, in step 912, the transfer function of the AR filters for auto power 
spectra of X s (t) and Yj(t) are computed as the square root of the prediction error (a z ) divided by the Fourier transform 
of the coefficients, i.e., 



15 



H Pz (f) = — — / z = X, Y. 



20 



i + X ci2 e " 

i = 1 



[0085] The auto/cross power spectral density values (PD(f)) is the magnitude of the complex product of H PX (f) and 
the complex conjugate of H Py </), i.e., 

25 

PC(f)= H px (f)*H py (f) 



30 PD(f) = IPC(f)l 

[0086] If requested, the same normalization used in step 808 is may be used here (on H P £f)). 
[0087] In step 914, the system checks whether the computation to be performed is a bispectrum computation, and 
if it is not, the system returns the requested auto/cross power spectral density array to monitor module 402. 
35 [0088] In steps 91 6, 91 8, and 920, the system sets the symmetries in the same manner as described above in steps 
814, 816, and 818. 

[0089] The system estimates the auto/cross bispectrum in steps 922, 924, and 926. The estimation process includes 
two major stages: the order selection and bispectrum computation. In step 922, two sequences of third-order moments, 
{R 3X (x)} and {R 3 y(x)} are computed using the following equation. 



H S 2 

RlzW = X X «(C)z3(frK), * = X, Y, and X = -L, 

45 1=1 t=Sl 



50 



55 



where s 1 = max (1 ,1-t), s 2 = min (N, N-x), and L is greater than the largest possible AR filter order (e.g. 50). 
[0090] In step 924, two matrices T x and T Y are formed as follows. 



Tt = 



(RlzirL) Riz{~L+l) ... K 3 z(0) ^ 
R3z(rL-l) Rlz(~L) ... K 3 z(-1) 

*3z(-2L) K 3 z(-2I+1) R3z(-L) 



z = X, Y 
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[0091] From the assumption we made about the AR filter of bispectrum, the orders and Oy of the AR filters of 
bispectra of X|(t) and Yj(t) are the ranks of the super matrices T x and T Y . Therefore, Ox and Oy are chosen by using 
singular value decomposition. Having found the orders, we obtain the coefficients of the AR filters of the bispectra by 
solving the following linear system of equations: 

5 

'*3z(0) *3*(1) *3z(O r ) > 

*3z(~D *3*(0) — *3z<Or-l) 

w 4 ' ... • 

Rlz(-Oz) R3zi-Oz+l) ... «3i(0) 
V J 

where the skewness (p z ) and the coefficients (b 1z , .... b 0zZ ), z = X or Y, can be obtained by solving the linear system 
15 of equations. 

[0092] The auto/cross bispectrum of X|(t) and Y|(t) are computed in step 926 as the cubic root of the triple product 
of the skewnesses (PxPyPy) 173 . divided by the triple product of the Fourier transforms of the AR filter coefficients (H z 
(0). i.e., 

BC(i v f 2 ) = tf x P v P y ) 1 °/H x (f 1 ) w x(W / M'i + y 



H Z (f ) = 1 + I »iz z = X,Y. 

i = 1 

30 and BR(f 1 ,f 2 ) is the real triple product for that same signal: 

BR(f l ,f 2 ) = P x (f 1 )*P x (f 2 )*P Y (f 1 +f 2 ) 

35 where the auto power spectra of Xj(t) and Y^t), P x (f) and P y (f), are computed by squaring the magnitudes of transfer 
function of the AR filters for auto power spectra of X^t) and Yj(t) (H PX (/) and H py </) ; respectively. If requested, the 
same normalization used in step 808 may be used here. Similarly, (P 2 ) 1/3 /f7 z (/) is divided by the square root of the sum 
of the square of its magnitude for certain frequency band, its largest peak value, or some similarly derived normalizing 
constant. 

40 [0093] After obtaining power spectrum and auto/cross bispectrum, the system computes the bispectral density array, 
the biphase, the bicoherence, and the square-rooted average real triple product (RTP) array in step 928 in the same 
way as in steps 826, 828, 830, and 832. In step 930, the system returns to the monitor module 402 the requested auto/ 
cross power spectral density array, bispectra! density, square-rooted real triple product, biphase, and bicoherence 
arrays. 

45 [0094] The above parametric equations used to compute the auto/cross bispectral arrays can be generalized to 
compute auto/cross higher-order spectral arrays. This will allow the computation of the trispectrum, quadspectrum, 
etc. Assuming that the arrays for a Kth-order spectrum are to be computed the following equations can be used: 
[0095] The auto/cross Kth-order spectrum: 

50 KC(f 1( f 2 f K . 1 ) = 

(<Px>*V'*' H x<W 
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Oz 

Hz (f ) = 1 + £ b i* e " 



z = X,Y. 



i = 1 



[0096] 



The real Kth-order product: 



KR(f v f 2 . .... f K .,) = P x (f t ) * P x (* 2 ) * - * P x( f K-i) * P Y ( f i +f 2 + - +f K-i) 
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[0097] After obtaining the auto/cross Kth-order spectrum, the system computes the auto/cross Kth-order spectral 
density array, the auto/cross Kth-order phase, and the auto/cross Kth-order coherence the same way as in the fre- 
quency-domain-based method. 

[0098] For illustration purposes Figs. 10(a) - 10(c) are graphs of sample autobispectral arrays showing frequency 
pairs 0 < f 1 < 30 Hz, and 0 < f 2 < 15 Hz. A bispectra! density array is shown in Fig. 10(a) where the Z axis represents 
the magnitude in decibels (db) of the coupled interaction between all appropriate frequency pairs f n and f 2 . Recall that 
the frequency pair (f v f 2 ) must adhere to the equation: 



where N = 256 Hz in this case. A bicoherence array is shown in Fig. 10(c) where the Z axis represents the normalized 
magnitude in percent (%) of the coupled interaction between all appropriate frequency pairs f t and f 2 . A biphase array 
is shown in Fig. 1 0(b) where the Z axis represents the phase in radians of the coupled interaction between all appropriate 
frequency pairs f, and f 2 . An array of square root of real triple product is shown in Fig. 10(d) where the Z axis represents 
the magnitude in decibels (db) of the coupled interaction between all appropriate frequency pairs f n and f 2 . 
[0099] Referring to figure 11 , a more detailed illustration of the diagnostic index generation module 410 will now be 
provided. In step 1 1 02, the system identifies the type of diagnostic assessment in progress. In a preferred embodiment, 
the 5 possible options are: 

1 . Depth of anesthesia, consciousness, pain & surgical stress. 

2. Cerebral ischemia and hypoxia. 

3. Cerebral intoxication (alcohol, narcotics). 

4. Evoked potential evaluation 

5. Cognitive process evaluation 

[01 00] In step 1 1 04, the system gets the auto/cross power spectrum, bispectrum, and/or higher-order spectrum arrays 
that are required for the computation of the requested diagnostic index using the sorting method described below. The 
various arrays that can be used in the generation of the diagnostic index are: auto/cross power spectrum; auto/cross 
bispectral density; auto/cross bicoherence; auto/cross bispectral real product; auto/cross biphase; auto/cross Kth-order 
spectral density; auto/cross Kth-order coherence; auto/cross Kth-order spectral real product; and auto/cross Kth-order 
phase; 

[0101] The sorting of auto/cross power spectrum, bispectrum, and higher-order spectrum arrays is an important 
feature as it provides a mechanism to compensate for changes in the energy distribution in these (and any other) 
spectra. The following is a general description of how the feature is supplemented: 

[0102] Based on an FFT derived from 4-second records, as described in the preferred embodiment, 120 data points 
can be computed for a power spectrum array that covers the frequency band 0-30 Hz (with 4-second records and a 
sampling rate of 1 28 samples per second, the resolution of the FFT is 0.25 Hz, and the range used is 30 Hz wide, thus 
there are 1 20 = 30 Hz / 0.25 Hz data points). When the 1 20 data points are sorted in descending order, the first element 
in the sorted array will correspond to the largest power spectrum value, and the last element will correspond to the 
smallest power spectrum value. A distribution histogram of the power can then be generated using the sorted array. 
The X axis on the histogram will represent power in dBs and the Y axis will represent the number of points in the sorted 
array that correspond to a particular X axis power value. If all points in the sorted array are added together, the sum 
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will represent the total power in the 0-30 Hz spectrum. If a number of adjacent points in the sorted array are added 
together, a portion of the histogram representing a percentage of total power is obtained. For example, in a particular 
EEG signal, the top 2 points in the sorted array represents the top 10% of the total power in the power distribution 
histogram. Similarly adding the bottom 70 points (for the same signal) in the sorted array will give the bottom 10% of 

s the total power in the histogram. Given this approach any portion of the power distribution histogram can be obtained 
by adding adjacent elements in the sorted array (top 25% of total power, middle 50% of total power, etc.) (given that 
one has empirically determined the transfer function from specific points to percentage of total power). 
[0103] By sorting, we are able to track regions of high activity and low activity (peaks and valleys) in the 0-30 Hz 
power spectrum without having to analyze specific narrow frequency bands. This is equivalent to mapping the power 

10 spectrum to its power distribution function and operating on fixed bands within that distribution function. This transfor- 
mation addresses some of the inconsistencies in the behavior of EEG power observed when hypnotic anesthetic agents 
are administered. More generally, the sorting scheme outlined above will transform any auto/cross power spectrum, 
bispectrum, and higher-order spectrum array of any dimension, into a one-dimensional distribution function of the 
values it contains. The one dimensional distribution is then divided into fixed bands that can be combined to produce 

15 a univariate diagnostic index. 

[0104] In step 1106, the reference auto/cross power spectrum, bispectrum, and higher-order spectrum arrays are 
sorted. The corresponding dependent arrays are re-ordered according to the sorted sequence of the reference array. 
A reference array is an array whose values are used as the primary sort key for a group of corresponding arrays that 
have the same number of variables and are identical in size to the reference array. For example, if the reference array 

20 were to have four elements and they were given the indices 1 , 2, 3, 4 before sorting, and after the sort the new order 
of the indices were 2, 1 , 4, 3 then one could use the same rearrangement to reorder any other array of the same size 
(in this case by placing the second element first, the first element second, etc.). In this way, one can use the sort of 
the reference array to rearrange the dependent arrays. In the preferred embodiment, the reference array is autobis- 
pectra! density and the dependent arrays are autobicoherence and the square rooted average real triple product. 

25 Autobispectral density was selected as the reference array because it provides information about the residual power 
at each frequency pair after random phase cancellations. Thus, the sort of the autobispectral density array provides a 
more stable means to select autobicoherence and real triple product values than would the sort of those arrays them- 
selves. A different array may be selected to satisfy other requirements. 

[0105] In step 1 1 08, the sorted auto/cross power spectrum, bispectrum, and higher-order spectrum arrays are each 
30 divided into bins as described earlier. The sum of the points in each bin for each array is computed and stored in a 
temporary variable. In step 1110, the clinically predetermined coefficient array for the desired diagnostic index is re- 
trieved from resident memory (or from the storage device). Each coefficient in the predetermined coefficient array 
corresponds to one of the temporary variables generated in step 11 08. In step 1112, the univariate diagnostic index is 
produced from the sum of all variables multiplied by their corresponding coefficients in the predetermined coefficient 
35 array. In step 1114, the program returns to the monitor module 402. 

[0106] The predetermined clinical coefficient arrays referred to above are essential to the device's ability to achieve 
clinically relevant diagnostic efficacy. The process adopted for generating these clinical reference arrays will now be 
described. Since a large number of possible reference arrays must be generated to accommodate all the diagnostic 
modalities of the system, only one will be discussed in detail. All other reference arrays are generated in a similar 
to fashion. For illustration purposes a method for generating the coefficients required to track depth of anesthesia using 
the derived signals FL and FR (of the preferred embodiment) is described below: 

[0107] In order to determine the clinical coefficients for a particular diagnostic index, raw data as well as clinical 
diagnoses are required. In the particular case described below, in order to develop an index which indicates anesthetic 
depth, EEG signals and assessments of the patient's response to clinical stimuli were collected. In one case below, 

45 the assessment is based on the patient's change in arterial blood pressure. In the other case, the assessment is the 
surgeon's judgment as to whether the patient had a motor-reflexive response. Once the data have been obtained, the 
various spectra and variables as described above may be generated. By performing statistical regressions on the 
processed data in conjunction with the clinical diagnoses, the coefficients which produce the index with the best pre- 
dictive diagnostic ability can be determined. 

50 [0108] In two separate studies EEG potentials were continuously recorded from a group of patients undergoing 
elective surgery. The recording period started at approximately 5 minutes prior to induction and lasted for the duration 
of the surgery. The derived signals FL, FR, PL, PR, FPL and FPR were acquired using the procedure described above. 
[0109] The purpose of the first study was to determine whether autobispectrum variables provide information about 
anesthetic depth at incision. Forty adult patients were studied. Anesthesia was induced with thiopental (up to 5.0 mg/ 

55 kg) and intubation performed after the administration of succinyl-choline. Patients were randomly assigned to receive 
isoflurane 0.75 MAC (Mean Alveolar Concentration), 1.00 MAC, or 1.25 MAC in 100% oxygen. End-tidal agent con- 
centration was monitored and after a period of steady-state had been achieved, purposeful movement in response to 
skin incision was assessed. Each patient was classified as either a "mover" or a "non-mover" based on the patient's 
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response to incision. 

[0110] The purpose of the second study was to determine whether autobispectrum variables provide information 
about predicting hemodynamic responses to laryngoscopy during induction with sufentanil or alfentanil. Forty adult 
patients were studied. Patients received premedication with oral diazepam (0.05 - 0.15 mg/kg) and were induced with 

5 thiopental (4.0 - 6.0 mg/kg) and 60% nitrous oxide in oxygen, followed by vecuronium (0.1 mg/kg) . Each patient was 
then randomly assigned to receive one of five dose regimens: normal saline; alfentanil 15 mcg/kg or 30 mcg/kg; suf- 
entanil 0.5 mcg/kg or 1 .5 mcg/kg. Laryngoscopy was performed 3 minutes after drug administration. Brachial blood 
pressure was measured every minute with a cuff device. Patients who exhibited a change in mean arterial pressure 
of more than 20% in response to intubation were classified as "responders"; those who did not exhibit such a change 

10 at intubation were classified as "non -responders." 

[0111] An autobispectral density, an auto bicoherence, and an auto square-rooted average real triple product array 
were generated for the derived signals FL and FR for each patient using a two minute period prior to the stimulus. The 
frequency band for which the bispectral arrays were computed was 0.25 - 30 Hz. Each bispectral array contained 3600 
data points. 

15 [0112] The resultant auto bispectral density, auto bicoherence, and auto square-rooted average real triple product 
arrays were sorted using the auto bispectral density array as the sorting reference array. The sorting was done using 
the algorithm described above. 

[0113] Eleven variables were produced from each of the sorted arrays as described below: 

20 Var1 =Sum of the largest 15 points in sorted array 

Var2 =Sum of points ranked 16th to 30th in sorted array 

Var3 =Sum of points ranked 31st to 50th in sorted array 

Var4 =Sum of points ranked 51th to 100th in sorted array 

Var5 =Sum of points ranked 101th to 150th in sorted array 
25 Var6 =Sum of points ranked 151th to 300th in sorted array 

Var7 =Sum of points ranked 301th to 500th in sorted array 

Var8 =Sum of points ranked 501st to 900th in sorted array 

Var9 =Sum of points ranked 901st to 1500th in sorted array 

Var10=Sum of points ranked 1501st to 2400th in sorted array 
30 Var1 1 =Sum of points ranked 2401 st to 3600th in sorted array 

[0114] The values of the 11 variables for each array were computed. As a result, there were 33 temporary variables 
per patient per signal. 

[0115] The 80 patients were then classified into two groups. The first group contained all the patients from the first 
35 study that moved at incision and all the patients from the second study that had a change in blood pressure of greater 
than 20% in response to intubation. The second had all patients from the first study who did not move at incision and 
ail the patients from the second study that had a blood pressure response of less than 20% for intubation. 
[0116] In order to produce a set of coefficients that would yield the most effective diagnostic index, a discriminant 
analysis was performed. The diagnostic index (l(c 0 , c n , - . c 33 )) for a set of coefficients (c 0 , c v . . c 33 ) is given by: 

40 

/(c 0 , c, c^)=c 0 HBlS A *c,+...+BIS K *c u ) 

+(S/C^*c 12 +...+B/C /c *c 22 )+(PS A *c 23 +...+PS /< *C3 3 ) 

45 

where BIS A through BIS K are the 11 sorted temporary variables from the bispectrum array; BIC A through BIC K are the 
variables from the bicoherence array; and PS A through PS K are the variables from the sorted square-rooted average 
real triple product array. The discriminant analysis, given the values of the temporary variables mentioned above and 
the responder/non-responder classification for each patient, produces the set of coefficients which yield the best sep- 
50 aration of responders and non-responders by the function I. Discriminant analysis algorithms are publically available; 
in this case, the ones used are from the statistics library available from IMSL (Houston, Texas) . Below is a sample list 
of coefficients generated using a database of 170 patients: 





for derived signals FL, FR 


Ci 


-4.28 
-0.65 
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(continued) 



5 



10 



20 





for derived signals FL, FR 




+0.57 


c 3 


+1 .21 


c 4 


-1.23 


c 5 


+2.63 


c 6 


-3.34 


c 7 


+2.11 


c 8 


+2.74 


c 9 


-3.08 


C 10 


0.0 


c„ 


-0.66 


C 12 


+0.04 


C 13 


-1.86 


C,4 


+0.50 


Cis 


+0.14 


^16 


-0.30 


C17 


+0.15 


C 18 


-0.08 


C 19 


-0.11 


C 20 


+0.05 


c 21 


+0.05 


C22 


-0.02 


C 23 


+0.67 


c 24 


-1.02 


C 25 


0.0 


°26 


.n iq 


c 27 


-1.27 


C 28 


+1 .20 


C 29 


+1.25 


C 30 


-2.15 


C 31 


-2.43 


C 32 


+3.16 


C 33 


+0.64 



[0117] For the two studies discussed above the univariate index was used to predict the response to the stimulus 
for each patient. The following is a summary of the results achieved: 

Sensitivity; predicting movement at incision = 96% 
Specificity; predicting no movement at incision = 63% 
Overall accuracy; predicting move/no move at incision = 83% 



Sensitivity; predicting >20% BP change at intubation =100% 
Specificity; predicting <20% BP change at intubation = 50% 
Overall accuracy; predicting BP change at intubation = 85% 

[0118] The example above shows one approach to obtaining a set of coefficients for a diagnostic application in a 
retrospective manner. Several other approaches can be used to separate the clinical populations being studied using 
a univariate index. Such approaches include but are not limited to linear regression, stepwise linear regression, logistic 
regression, and stepwise logistic regression. Of course regardless of which method is used to retrospectively compute 
the coefficients, performance of the final index must be confirmed in a prospective trial prior to using it in patient care. 
[0119] The analytic process described above is used to generate the reference databases for cerebral ischemia, 
cerebral hypoxia, consciousness, degrees of intoxication, altered evoked potential responses, and normal or abnormal 
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cognitive processes including but not limited to identifying patients with Alzheimer's disease and HIV-related dementias. 
[0120] In addition to quantifying the depth, adequacy of anesthesia, and pain responses during surgical stress, the 
system and method of the present invention may also be used to assess a myriad of cerebral phenomena that alter 
the nonlinear frequency structure of the EEG as quantified by bispectrum and higher-order spectrum approaches. Such 
cerebral phenomena include but are not limited to, cerebral ischemia, cerebral hypoxia, level of consciousness, degree 
of cerebral intoxication, altered evoked potential responses, and normal or abnormal cognitive processes caused by 
neurological disorders like Alzheimer's disease or HIV-related dementias. 

[0121] Although power spectrum and bispectrum analysis techniques have been applied to the EEG signal for diag- 
nostic purposes, as was discussed in the background above, higher-order spectral approaches have never been used. 
Furthermore no power spectrum, bispectrum, or higher-order spectrum technique has ever been used in conjunction 
with the sorting method described above. Specifically, the system and method of the present invention sort various 
auto/cross power spectrum, bispectrum, and higher-order spectrum arrays, divides the sorted arrays into bins, sums 
the variables in each bin, multiplies the value from each bin by a clinically derived coefficient, and finally adds all 
variables together to generate a univariate diagnostic index. The different arrays that can be used are: auto/cross 
power spectrum, auto/cross bispectral density, auto/cross bicoherence, auto/cross biphase, auto/cross average real 
triple product, auto/cross Kth-order spectral density, auto/cross Kth-order coherence, auto/cross phase, and auto/cross 
real product. 

[0122] While the foregoing invention has been described with reference to its preferred embodiments, various alter- 
ations and modifications will occur to those skilled in the art. All such alterations and modifications are intended to fall 
within the scope of the appended claims. 

Claims 

1. A method of generating a diagnostic index for quantifying the presence or absence of a phenomena which is 
represented by spectral values generated from acquired electrical signals of the phenomena, said method com- 
prising the steps of: 

sorting said spectral values into predetermined bins of ranges of spectral values; 
summing all of the spectral values in each bin; 

multiplying the sum of spectral values in each bin by a predetermined coefficient to obtain a bin product; 
summing said bin products to obtain a diagnostic index which represents a degree of presence or absence of 
said phenomena. 

2. A system for generating a diagnostic index for quantifying the presence or absence of a phenomena which is 
represented by spectral values generated from acquired electrical signals of the phenomena, said system com- 
prising: 

means for sorting said spectral values into predetermined bins of ranges of spectral values; 
means for summing all of the spectral values in each bin; 

means for multiplying the sum of spectra! values in each bin by a predetermined coefficient to obtain a bin 
product; 

means for summing said bin products to obtain a diagnostic index which represents a degree of presence or 
absence of said phenomena. 
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Fig. 9 
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